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I. INTRODUCTION 

Large fluctuations play an important role in many 
physical phenomena, an example being spontaneous 
switching between coexisting stable states of a system, 
like switching between the magnetization states in mag- 
nets, or voltage/current states in Josephson junctions, 
or macromolecule configurations or populations. Typi- 
cally, large fluctuations are rare events on the dynamical 
time scale of the system. A theoretical analysis of such 
events goes back to Kramers 1], who considered the rate 
of switching of a Brownian particle from a potential well. 
The problem of the switching rate and the probability 
distribution becomes more complicated for systems away 
from thermal equilibrium, as in this case the probability 
distribution is no longer of the Boltzmann form. A rig- 
orous mathematical approach to the problem was devel- 
oped and many results have been obtained for dynamical 
systems without delay driven by white Gaussian noise 
and for Markovian reaction and population systems, cf. 
[MTU]. More recently the problem of large rare fluctu- 
ations in white-noise driven systems with delay was ad- 
dressed in the mathematical literature |I1| . 

Delay naturally arises in dissipative dynamical sys- 
tems. In such systems, dissipation results from the cou- 
pling to a reservoir: motion of the system causes changes 
in the reservoir, which in turn affect the motion. The un- 
derlying reaction of the reservoir is generically delayed. 
Along with the dissipative force, the reservoir exerts a 
random force on the system. If dissipation is delayed, 
the random force has a finite correlation time. These ef- 
fects have been attracting much attention in the context 
of optomechanics and the dynamical back-action ^12:. 

In this paper we develop a formalism for studying large 
rare fluctuations in classical systems with delayed dissi- 
pation. Much of the analysis refers to the case of Gaus- 
sian noise, but the results can be immediately extended 
to non-Gaussian noise as well. An important part of 
the paper, that allows us to test the general formula- 
tion, is the analysis of coupling to a reservoir in thermal 



equilibrium, where the noise and the dissipative force are 
connected by the fluctuation-dissipation relation |I3j . 

Central for the analysis is the idea of the optimal fluc- 
tuation. In a large rare fluctuation the system is brought 
from its stable state to a remote state in phase space. 
This requires a large deviation of the noise from its root 
mean square value. Different noise realizations can result 
in the same outcome, but they have different probability 
densities; for Gaussian noise, the difference is exponen- 
tial. The overall probability of a large fluctuation of the 
system is determined by the most probable, or optimal 
appropriate realization of the noise. 

As a consequence, in a fluctuation to a remote point 
in phase space or in switching the system is most likely 
to move along a well-defined (optimal) trajectory, that 
corresponds to the optimal noise realization. Using the 
approach [14' and its extensions, the narrow peak in the 
distribution of the trajectories has indeed been seen in 
simulations and in the experiments, cf. [T5U20| . 

An important feature of large rare fluctuations in 
Markovian (no delay) systems in thermal equilibrium 
with a bath is that the optimal fluctuational path is the 
time reversed path in the absence of noise, cf. [ST . This 
can be understood from the argument that, in relaxation 
in the absence of noise, the energy of the system goes 
into the entropy of the thermal reservoir, whereas in a 
large fluctuation the entropy of the reservoir goes into 
the system energy. The minimal entropy change corre- 
sponds to a time-reversed process [T3]. In other words, 
the optimal trajectory for a large fluctuation corresponds 
to the noise-free trajectory with the inversed sign of the 
friction coefficient. One can view this property also as a 
consequence of the symmetry of transition rates in sys- 
tems with detailed balance discussed for diffusive systems 
described by the Fokker-Planck equation by Kolmogorov 
[22] (optimal fluctuational paths and the path distribu- 
tion were not discussed in Ref . . 

We show below that the situation is more complicated 
if the dissipative force is delayed. We consider linear 
coupling to a thermal reservoir, which leads to a delayed 
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viscous friction. The model is described in Sec. ITTl below. 
A variational problem for finding optimal fluctuational 
paths in systems with delayed relaxation in the presence 
of Gaussian noise is formulated in Sec. |III| Both the prob- 
lems of the tail of the probability distribution and escape 
from a metastable state are considered. In Sec. IIVI we 
show that, if the noise has thermal origin, the tail of the 
probability distribution remains to be of the Boltzmann 
form in the presence of dissipation delay. We also con- 
sider the time-reversal symmetry of the most probable 
trajectories. In Sec. |V] the results are illustrated using 
an exponentially correlated thermal noise. In Sec. |VI| we 
give explicit expressions for the logarithm of the prob- 
ability distribution and the escape activation energy for 
the case where the correlation time of the noise and the 



dissipation delay are short. Sec. VII contains concluding 
remarks. 



II. A SYSTEM LINEARLY COUPLED TO A 
THERMAL BATH 

We start by sketching the derivation of the equation 
of motion with delayed dissipation for a particle coupled 
to a thermal bath. This problem was addressed back 
in the mid-60s [SS", [24l. In contrast to Refs. 23' and "241, 
our analysis refers to a classical particle in a generally 
nonparabolic potential U{q), where q is the particle co- 
ordinate. The particle has a unit mass and the bath has 
a quasi-continuous frequency spectrum. The coupling is 
linear in the particle coordinate q. The total Hamiltonian 
of the system and the bath is 



(1) 



Here, p is the momentum of the particle; i?b is the Hamil- 
tonian of the bath in the absence of the interaction; /ib 
is a function of the dynamical variables of the bath only. 
A familiar example is a bath that consists of harmonic 
oscillators, with h^, linear in the oscillator coordinates qk, 
cf. 
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hh = ^£kqk- (2) 



The analysis below is not limited to this model; it imme- 
diately extends also to the case where the interaction Hi 
is linear in the particle momentum p. 
In the equation of motion of the system 



m + U' {qit)) + h^it) = 0. 



(3) 



function h\^{t) itself depends on q{t') with t' < t, because 
the bath is perturbed by the system. In much of the 
paper we assume the interaction to be weak, i.e., that 
is proportional to a small constant. Then the response 
of the bath to the motion of the system can be described 
using the generalized susceptibility ah- It determines 



the mean value of hh if, instead of the coupling to the 
considered dynamical system, the bath were driven by a 
time-dependent force F{t), with energy —F{t)hb. In the 
considered case the role of F{t) is played by —q{t), 



dt'ah{t-t')q{t'), 



(4) 



where we assumed that the interaction was adiabati- 
cally turned on at t — )■ — oo. In the model ^ ah{t) = 
J^ki^k/^k) sin ujkt [25l [28]; Eq. Q applies in this case 
for an arbitrarily strong coupling. 

In Eq. h^\t) is the value of h\^(t) in the absence 
of the interaction with the system. It is a random func- 
tion of time; for example, in the model ^ the random- 
ness comes from the randomness of the amplitudes and 
phases of the non- interacting oscillators [5S]. We set 
{h^^\t)) — 0. The power spectrum of h^^"^ can be written 
as 2Re $/i(a;), where 



dtexp{iiot)(j)fi{t), 



(5) 



The power spectrum is related to the susceptibility ah 
by the fluctuation-dissipation theorem [13]. 

If hi, is nonlinear in the dynamical variables of exci- 
tations of the bath, Eq. Q can be obtained by calcu- 
lating the first-order in Hi correction to these variables 
and expanding hb to the first order in this correction. 
Function ah{t — t') is determined by the mean value of 
the coefficient at q{t') calculated in the absence of the 
system-to-bath coupling. The remaining random part is 
of higher order in the interaction than h'^^ and therefore 
is disregarded. As will be shown separately, the decou- 
pling that leads to Eq. Q applies if the correlation time 
of fluctuations of the bath is small compared to the effec- 
tive relaxation time of the system, which is determined 
by the interaction Hi. For the coupling ^ there is no 
this limitation. 

Using that, from the fluctuation-dissipation theorem, 
ah{t) = -P(l)h{t) for t > (/3 = l/ksT), one obtains 
from Eqs. ([l]) and Q the equation of motion of the sys- 
tem coupled to the bath in the form 



q{t) + U'h{q{t)) + P / dt'cj,h{t')q{t - t') = f{t), 
Jo 

Uh{q) = U{q) - (/?/2)g20^(O). 



(6) 



Here f{t) is a random force; if the only source of this force 
is the coupling to the bath, then f{t) = h^^\t). However, 
often in the experiment noise comes from external sources 
that are not in thermal equilibrium with the system, and 
f{t) accounts for such noise. 

The integral term in Eq. ^ describes dissipation of 
the system due to the coupling to a thermal bath. The 
dissipation is delayed. Within the linear response ap- 
proximation it is linear in the velocity of the system, but 
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depends on the velocity history. The coupUng to the 
bath leads also to the renormalization of the potential 
of the system U{q) — >■ Uh{q)- We note that it is natural 
to count q off from a minimum of the "bare" potential 
U{q) (a constant shift of q can be incorporated into H^,). 
Then the renormalization (|6| corresponds to softening of 
the potential near this minimum, since 4>h{Q) > 0. 

A. Stationary states 

The stationary states of the system ([6| in the absence 
of noise are located at the extrema of Uh{q)- Near its 
extremum go the potential Uh can be linearized in 6q = 
q — qo and the solution of Eq. Q can be sought in the 
form 6q{t) oc exp(Ai) with A given by equation 

A2 + t/;'(Qo)+^/3A$,,(zA) =.0. (7) 

Equation ([T]) is simplified if the coupling to the bath is 
so weak that the last term is a perturbation compared to 
U'fl^. In this case, if Uh has a minimum at qo = qa, with 
U'^iqa) ^ a;2 > 0, 

Xa± « ±^Wa - r, r = (7r/3/2)Re <^h{oJa), (8) 

where = + (7r/3/2)Im $;i(a;a). Equation ^ applies 
provided /3|<I>ft(wa)| < Wq. 

Since the thermal noise power spectrum 2Re ^h{^) is 
non-negative, F > 0. Therefore in the presence of delay 
a minimum of Uh still corresponds to an asymptotically 
stable state of the system, an attractor. Parameter F is 
the characteristic relaxation rate of the system. Since 
r ^ Wa, the stable state is a focus on the phase plane 
{q,p) and the motion near {qa,Pa = 0) is underdamped. 

Generally, systems with delay have an infinite- 
dimensional phase space. Therefore Eq. ([7| should have 
more than two solutions. However, function is an- 

alytical for Im u! > 0. Then <i>;i(iA) does not diverge for 
Re A > 0, and if the coupling is small, the last term in 
Eq. remains small for Re A > 0, which means that 
there are no solutions of Eq. Q with Re A > in the 
weak-coupling limit. In turn, this means that a minimum 
of Uh{q) is an attractor. 

In many cases of physical interest the correlator (phit) 
exponentially decays for large times, 4>h{t) on exp{—t/tc) 
for t -4- oo, where tc is the correlation time of bath fluc- 
tuations. Then for a very weak coupling Eq. ^ has 
a root Re A « —1/tc, which describes a comparatively 
fast relaxation of the bath when the system is close to 
the attractor, Ttc ^ 1. We will not discuss in this pa- 
per the case of a power-law decay of correlations in the 
bath, which has attracted much attention in the context 
of quantum tunneling 29 , although some of the results 
apply to this case, too. 

A local maximum of Uh{q) is a saddle point; we will 
use the notation qg for this point. Within a naive pertur- 
bation theory, near the saddle point in the limit of small 
the eigenvalues are 

A5± « ±ns - {^p/2)^h {±^s)) , (9) 



where Q.s = \U"{qs)\^^'^ ■ The coupling-induced renor- 
malization of the root A5+, that describes moving away 
from qs, is small for small \^h\- However, the change 
of the root A^^ can be significant even in the small- |<I>/i| 
limit. Indeed, if (j)h{t) oc exp(—t/tc) for i — > 00 and 
t-^ < \U'^{qs)\^/^, the correction oc <^h to A^- in Eq. ^ 
diverges: formally, ^h(j-^) diverges for — Im lo > 1/ic- 
Physically, the system moves too quickly for the bath to 
follow it. In this case Re A^- approaches —1/tc in the 
weak-coupling limit, i.e., an arbitrarily weak coupling to 
the bath leads to a finite change of the "stable" eigen- 
value near the saddle point. 

A dramatic change of the dynamics can occur also for 
the coupling of the form ([2]), if the bath frequencies form 
a finite- width band [30^. We will not consider this case 
here. 

Even where the coupling to the bath is weak in the 
sense that the decay of the system is slow compared to the 
decay of correlations in the bath, the relaxation rate of 
the system can exceed the frequency Wq. In this case, for 
^hi^) smooth near io — 0, the motion near the potential 
minimum is overdamped, 

Xa+~-ijl/nl3<Ph{0), A,_ « -7r/3$,,(0), (10) 

with |Aa-| ^ uja ^ \^a+\- Obviously, the eigenvalues 
Xa± are real and negative, indicating that the potential 
minimum remains an attractor. We call the correspond- 
ing parameter range the small inertia limit. Indeed, if 
we incorporate the particle mass m into the Hamiltonian 
Hq and define U'j[{qa) — muj'^, we see that Xa- oc 1/to, 
whereas Xa+ oc to. The root Xa+ characterizes the slow 
motion of the particle coupled to the bath, other degrees 
of freedom relax much faster. 

Similarly, near qs for small inertia we have 

A5- « -^/3$h(0), Xs+-\U'^:iqs)f/7rP^h{0) (11) 

with |A5_| > \U^{qs)\ > IA5+I. The potential maxi- 
mum remains a saddle point, with motion away from it 
being slow compared to the motion toward it. Approach- 
ing the saddle point is characterized, essentially, by the 
relaxation rate of the bath when the system is at the sad- 
dle point. Of relevance to the motion of the system in 
the limit of small mass is primarily the root Xs+ ■ 

If the coupling is described by the model ([2| and is 
not weak, we have not found a simple explicit expression 
for the eigenvalues A. However, one can expect that the 
minima of Uh{q) remain stable states. This is a conse- 
quence of the condition Re <f>/i(a;) > for Im oj — 0. 
Indeed, because of this condition Eq. ([7| has no roots 
with Re A = for go = Qa- Therefore, since for weak 
coupling Re Aa < 0, and given that the dependence of Xa 
on the coupling strength is smooth, as Eq. Q suggests, 
as we increase the coupling strength the roots Aa will 
never cross the axis Re A^ = 0. Hence, the state (gajPa) 
will remain stable, see an example in Sec. [Vl 
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III. VARIATIONAL PROBLEM FOR OPTIMAL 
FLUCTUATIONS 

Large rare fluctuations in systems with delayed dis- 
sipation and with non-white noise can be analyzed by 
extending the approach developed for systems with no 
delay. The underlying idea is that, before a large fluctu- 
ation occurs, the system performs small-amplitude fluc- 
tuations about its initially occupied stable state. These 
small fluctuations persist for a long time that largely ex- 
ceeds the relaxation time. When ultimately there oc- 
curs a large fluctuation to a given point {q,p) in phase 
space, the system is most likely to move along the op- 
timal trajectory that corresponds to the most probable 
appropriate realization of the noise, as outlined in the 
Introduction. 

If the noise f{t) in Eq. ([6| is Gaussian and station- 
ary, with zero mean and correlation function (j)f[t~t') ~ 
{f{t)f{t')), the probability density functional of noise re- 
alizations is Vflfit)] = cxp(^-7^/[/(^)]/D) [27], where 



%[/(*)] = 4 / dtdt'f{t):F{t-t')f{t'), 



dtiFit - ti)(f>f{ti - t') = 2Dd{t - t'). (12) 



Here, D is the characteristic noise intensity; in the case 
of thermal noise, f{t) = h^\t) and D = kgT. We as- 
sume that D is small, so that on average the amplitude 



of fluctuations of the system about its attractor is small. 
Function F{t)/2D is the inverse of the noise correlator 

To logarithmic accuracy, the probability density of 
reaching a state {q^p) is 



p{q;P) — const X exp [—R{q,p)/D] , 



(13) 



where R{q,p) is the minimum of the functional TZf with 
respect to the noise trajectories that bring the system 
from the attractor to the state {q,p)- As for Gaussian- 
noise driven systems without delay [31], R{q,p) is given 
by a solution of the variational problem 

i?(g, p) = min [f{t)] + f dt x{t) [q{t) + {q{t)) 



dt'Mt-t')q{t')-f{t) 



(14) 



where the minimum is taken with respect to functions 
q{t), f{t), and x(i). The auxiliary function x(i) is a La- 
grange multiplier: it accounts for the interrelation be- 
tween the trajectory of the system q{t) and the force 
trajectory /(t), which is given by Eq. ^. 

Equation ( 14 ) describes the optimal fluctuational tra- 



jectories of the system and the noise. The boundary con- 
ditions for the trajectories where the system arrives from 
the attractor {qa,Pa = 0) occupied for t — >■ — oo to a given 
{q,p) at a given time t (we set t = 0) are 



fit) ^ 0, t 



fit),x{t) 

oo; xit) 



►0, q{t) 
0, i > 0; 



Qa, p{t) 
q{t - 0) 



'Pa 

q, 



for t 

p{t^ 



0) - p. 



(15) 



r 



The boundary condition for i = corresponds to the pic- 
ture of optimal fluctuation to a point (g,p) in which, once 
the system has reached this point, its further dynamics is 
no longer relevant. Respectively, the force should evolve 
for t > so as to minimize TZf independent of how its 
evolution affects the system. This is formally described 
by setting x = for t > 0. Alternatively, one can set the 



upper limit of the integral over t in Eq. ( 14 ) equal to zero, 
see below. Since 7?./ is positive definite, its minimum is 
reached for f{t) = 0, and therefore / — )■ for t — >■ — oo. 



The boundary condition (151 for t 



— oo allows for 
the fact that the system starts from the attractor and the 
optimal value of the noise is zero. A better understanding 
of this condition comes from the analysis of the equations 
of motion on the optimal trajectory. 



A. Optimal trajectory 



It follows from Eqs. ( 14 1 that on the optimal trajectory 



dt'T{t-t')!{t') = x(t), 



(16) 



x{t) + K{q{t))x{t) - p / dt'<t>h{t' - t)x{t') = 



while the interrelation between f{t) and q{t) is of the 
form ([6]). Equation ( 16 ) for x directly applies to the prob- 



lem of escape discussed in Section[III B[ In the problem of 



reaching a given state, as seen from Eq. (15), x(0 is dis- 
continuous at the instant t — when the system reaches 
the targeted state. Therefore x(^)j along with a smooth 
part at t < and t > 0, has a term —xi—£)S{t — e) with 
£ ^ +0. 



If in the variational problem ( 14 1 the integral over t 
(the one with the integrand oc x(0) is taken from — oo to 
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the instant of observation t = 0, Eq. (16) for x reads 



m + K{q{t))x{t) - P / dt'^^,{t' ~ t)x{t') 



t < 0. 



(17) 



It coincides with Eq. ( 16 ) for x is one aUows for the afore- 



mentioned (5-function in 



From Eqs. ( 12 ) and ( 16 I 



D 



dt'(^f{t-t')x{t'). 



(18) 



An interesting and somewhat counterintuitive feature 



of Eqs. ( 16 1 and pT| is that the time evolution of x{i) is 



acausal, the value of x(*i) depends on xi^'i) with t'l > ti. 
This does not make the equations iU-defined, since we 
are solving a boundary-value problem, where we know 
where the system arrives ai t ~ and what happens to 
the noise and x{t) after that. 

Equations ([6]) and ( 16 ) are simplified near the attrac- 
tor, where U'fXq) ~ w^(g — qa) and the equations become 
linear. In particular, for weak coupling to the bath x(t) 
for t — > ~oo has the form 



X(^) = X- exp{-Xa-t) + x+ cxp(-Aa+t), 



(19) 



where Aa± are given by Eqs. ([8|) or ( 10 ) and x± are ar- 
bitrary constants. Since Re Xa± < 0, the solution (19) 



satisfies the boundary condition (15) for t — > — oo. Gen- 



erally, there are also other terms in x(Oj see below, but 
they decay faster than ( 19 1 if the correlation time of the 



bath fluctuations is small compared to |Re Xa±\ 



assumed in Eq. (19). 

By comparing Eqs. (16) and ([7| one can see that, if 
A is a decrement that characterizes the approaching of 
the coordinate to the stable state in the absence of noise, 
then —A describes the increase of x{t) near the stable 
state for t — >• — oo. 

If the correlation time of the noise f{t) is shorter than 
the relaxation time of the system, then asymptotically 
for t — > — oo 



fit) 



fi,ex.p{-Xaut), 



*/(-*A„±)] 



(20) 



Here, ^fi^^) = tt^^ dt(l>f{t)exp{iujt); function 
2Re ^f{uj) is the power spectral density of the noise 
f{t). Equation (20) shows that on the optimal trajectory 
f{t) — > for t — oo, in agreement with the boundary 



condition ( 15 ). 



The deviation of the coordinate and momentum of the 
system from {qa,Pa) is also exponential in time for t — > 
— oo; from Eqs. ^ and (20) for a short noise correlation 
time 



lit) - 9a= - ^fv CXp(-Aaj,i) 



X {7r/3Aa. [<^h{iXa.) + '^h{-iXa.)]} ' . (21) 



From Eqs. m 



(21), the optimal trajectory is fully 
determined by the parameters x- and x+- These pa- 
rameters must be found from the boundary condition for 
t = 0. 



B. Escape problem 

Noise can also lead to escape of the system from the ini- 
tially occupied attractor and switching to another attrac- 
tor. To find the probability of the corresponding large 
fluctuation per unit time one should again minimize the 
functional TZf with respect to noise realizations that lead 
to escape. The key here is to note that, after the noise 
f{t) and the memory kernel ijihit) will have decayed, the 
system should be outside the basin of attraction of the 
initially occupied attractor or at least on the boundary 
of this basin. 

A correlated noise decays in time smoothly, it takes 
infinite time to decay to zero, as is also the case, generally, 
for the memory kernel. On the other hand, for t ^ oo the 
system approaches a stationary state. From Eq. ( 19 ) such 



state may not be an attractor, since x(i) would diverge 
there for i — >■ oo. Therefore it must be a saddle point g^, 
a local maximum of the potential Ut (?) ■ 

From the above arguments, which are similar to the 
corresponding arguments in systems without delay |31| . 
it follows that, to logarithmic accuracy, the escape rate 
has the form 



I^e(xexp(-i?^/D), 



(22) 



where the effective activation energy Ra is given by the 
solution of the variational problem ( 14 ) with the bound- 
ary condition for t — !■ — oo of the form of Eq. (15), 
whereas the other boundary conditions are 

/(t),x(i) -^0, q{t)^qs for t oo. (23) 



It follows from Eq. (16) linearized near qs that, for 
weak coupling, x(i) decays as exp(— A^+t) for t — >■ oo. 
If the correlation time of /(t) is smaller than I/A54., on 
the optimal escape trajectory f{t) and q{t) — qg decay 
in the same way. Otherwise the decay of f{t) and q{t) 
is controlled by the decay of f{t) as given by Eq. (18) 
or the decay of (phit). This completes the formulation 
of the problem of the activation energy of escape in the 
presence of delay. 

We note that the explicit expression for f{t) in terms of 
x(t), Eq. (18), allows eliminating /(t) and reducing the 



variational problem ( 14 ) to that for two coupled func- 
tions, x(^) and q{t). In fact, it corresponds to integrat- 
ing over realizations of f{t), and the resulting functional 
(besides the g-dependent part) is the characteristic func- 
tional of the noise 'P/[«x(0] [13 ■ The formulation can 
be extended also to non-Gaussian noise, and there, too, 
one can use the characteristic functional of the noise; 
the boundary conditions for the optimal fiuctuation to a 
given state of the system and for escape are the same as 
for Gaussian noise. 
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IV. LARGE FLUCTUATIONS INDUCED BY 
THERMAL NOISE 



Optimal trajectories of the system and the noise can be 
found and the logarithm of the probability distribution 
can be obtained in an explicit form in the case where 
the noise f{t) is of purely thermal origin, f{t) — h'^\t). 
In this case = </>/i(t) and the noise intensity D = 

UbT =1/ p. One can show that the variational equations 
for q{t), f{t), and x{t), Eqns. ^ and for the 
time prior to reaching the "target" point {q,p) at t = 
have a solution: 



t < 0, 



q(t) + U'„{q{t)) - f3 / dt'Mt' - t)q{t'); (24) 



function f(t) is given by Eq. (181 with x{t) — q{t) for 
t < and xit) = for i > 0. 



The solution ( 24 ) applies also to the problem of escape, 



except that the integral in the right-hand side of the sec- 
ond equation goes from i to oo and the solution applies 
in the whole range —oo < t < oo. 



Noting that in the variational problem (14) TZf[f{t)] — 
I I^oo ^tx{t)f{t) and expressing /(<) in this expression in 
terms of q from the equation of motion of the system ^ , 
one obtains 



Ra 



U,M - Uhiqa), 
Uhiqs) - Uh{qa)- 



(25) 



Thus, as expected, in thermal equilibrium the probabil- 
ity distribution of the system ( 13 1 is of the Boltzmann 



form, with a renormalized potential due to the interac- 
tion (however, our calculation does not give the prefactor 
of the distribution). The activation energy of escape is 
given by the renormalized height of the potential barrier. 
We now discuss the form of the optimal trajectories 



(24) for (5-correlated fluctuations of the thermal bath, 
where (puit) = ^TkBT5{t), with F being the coefhcient of 



viscous friction. In this case the trajectory q{t)^ Eq. (24), 



differs from the trajectory q{t) in the absence of noise, 
Eq. ([6| with f{t) = 0, by time inversion, t —t, q — >• —q, 
cf. [21 . Using the "prehistory problem" formulation 
this symmetry was seen in simulations and experiments. 
As mentioned in the Introduction, it can be understood 
by noting that in a fluctuation the system gets energy 
at the expense of the decrease of entropy of the ther- 
mal bath. Taking a time-reversed path minimizes the 
required entropy change. 

Care must be taken when applying the time-reversal 
argument in the general case of delayed dissipation. In- 
deed, to define a noise-free trajectory it is insufficient to 
specify the starting point in phase space {q,p); instead 
one has to specify the whole history of motion before the 
state (q,p) was reached. 



A. Symmetry of the most probable trajectories 

The previous analysis referred to the comparison of the 
most probable trajectories that come to a given state as 
a result of a fluctuation and the trajectories that start 
at this state and in the absence of noise go to the at- 
tractor. Instead one can consider trajectories that go 
through a given state. They reach it from the vicinity 
of the attractor as a result of a fluctuations, and then, 
as the fluctuation decays, they go back to the attrac- 
tor. The distribution of such trajectories was studied for 
Markovian systems by Luchinsky and McClintock [15j . 
It was shown that, for an overdamped system in ther- 
mal equilibrium, the distribution peaks at a trajectory 
that is symmetric with respect to time reversal, if time is 
counted off from the instant of reaching a chosen state. 
In other words, the segments of the most probable trajec- 
tory, which correspond approaching the state and moving 
away from it, are symmetrical. 

We now consider the most probable trajectory that 
goes through a given state {q,p) for systems with de- 
layed dissipation. We assume that the trajectory passes 
this state at t = (obviously, this instant is arbitrary). 
The segment of the trajectory for i < is the optimal 
fluctuational trajectory. In thermal equilibrium it is de- 
scribed by Eq. (24). The most probable trajectory after 
the state {q,p) has been reached is described by Eq. (|6|. 
In contrast to Markov systems, the force f{t) in this equa- 
tion is generally nonzero, because the noise is correlated. 
Since it was nonzero for t < 0, where it drove the system 
against the potential gradient, it does not instantly go 
to zero for t > 0. The most probable value of f{t) for 
t > is given by Eq. (flSl) with D = ksT. This value 



maximizes the probability density of a realization of f{t) 
that brings the system to state {q,p). 

From Eqs. ^ and ( 18 ) with account taken of the rela- 
tions x{t) = q for t < and x(i) = for i > 0, we obtain 
that the most probable trajectory of the system after the 
state {q,p) has been reached is described by equation 



-/? / dt'<j}h{t~t')q{t') [t > 0). (26) 
10 



A comparison of Eqs. ( 24 ) and ( 26 ) shows that, in ther- 



mal equilibrium, for the most probable trajectory that 
goes through a state {q,p) with p ~ 0, the segments that 
lead to the state from the attractor [Eq. (24)] and from 
the state to the attractor [Eq. (26)] have time-reversal 
symmetry t — >■ —t, q — > —q. This symmetry is illustrated 
in Fig. [l] below. 

We emphasize again that, in contrast to Markov sys- 
tems, one cannot compare separately the most probable 
trajectory to a given state and the most probable tra- 
jectory from a given state, because the latter is not de- 
fined. Only a single optimal trajectory that goes through 
a given state has the symmetry. This trajectory is fully 



determined by the state: the integrals in Eqs. (24) and 
(26) go to/from t = 0, where the system is at the given 
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{q,p)- Another distinction from Markov system is that 
the trajectory is smooth at {q,p) even where p ^ 0. 



V. THERMAL NOISE WITH EXPONENTIALLY 
CORRELATED POWER SPECTRUM 

To illustrate the general results we will discuss the case 
where the noise spectrum has a simple form, 



(t)hit) = Ch^cksT exp{-}t\t\), 



(27) 



where >f = 1/tc is the reciprocal correlation time and Ch 
characterizes the noise intensity; for the coupling model 
([2]) C is independent of temperature. The limit >f — 00 
corresponds to (S-correlated noise. 

For the model (27 1, Eq. ([t]) for the eigenvalues that 
characterize noise-free motion near attractor takes the 
form 



{X + >t){X^ + Lul) + Ch>iX = 0. 



(28) 



If the coupling is weak, two of the roots of this equation 
Ai,2 are given by Eqs. ^ for Chx'^ + i^a) ^ or 
Eq. (|lO| for X ^ Ch ^ i^a, respectively. In addition, 
Eq. ( |28[ ) has a root with a much larger (in the absolute 
value) real part 



A, 



(-A3 > -Re Xa±). 

In the opposite limit Cu 3> >f 3> w^, which can be 
of relevance for the model ([2]), the roots become Ai w 
-i^l/Ch and A2,3 ~ - [x T i^'^xCh - x'^] l'2.. For Ch ~ 
X ^ bJa the real parts of all 3 roots Ai^2,3 are of the same 
order of magnitude. It is easy to see that |Re A1.2.3I < >f, 
i.e., correlations of the bath decay faster than q{t). 

An explicit solution of Eq. (24) for the optimal tra- 



jectory to a state {q,p) can be obtained for a harmonic 
potential Uh{q) — w^g^/2. It reads 



q{t) =Y,(1^<^~^'\ Y.X,{x+X,)-\,^Q, 

i i 

^g, ^A,g, = (i- 1,2,3). (29) 

i i 

For the coupling to the bath of the form ([2| , this solution 
applies for an arbitrary coupling strength. It shows that 
the optimal fluctuational trajectory is a superposition of 
three exponentials. 

The most probable trajectories for reaching a given 
state {q,p) and then moving back to the attractor for 
the model (271 are shown in Fig. [I] The sections of the 
trajectories to and from the targeted state are symmetric 
for p = Q. The symmetry is lost if p 7^ 0, but the overall 
trajectory to and from the state is smooth. In contrast, 
for Markov systems, on the most probable trajectory that 
reaches a state {q,p) the derivative dp/dq is discontinuous 
at this state for p ^ 0. However, in such systems the most 



probable trajectory to state {q,p) has a symmetric noise- 
free counterpart that starts from (g, —p) and goes to the 
attractor. Such counterpart is not generally defined for 
systems with delay. 
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FIG. 1. The most probable trajectories of thermal fluc- 
tuations in which the system reaches targeted states (q,p) 
(the solid lines). The dashed lines show the sections of 
the trajectories after the targeted states have been reached. 
Plots (a) and (b) refer to a system in a parabolic poten- 
tial, Ur ~ and the exponentially correlated noise. 



Eq. (271, with x/u)a = l,Ch/wa — 1.5. The targeted states 
are (g = 2,p = 0) in (a) and (g = 2,p/u)a = 1) in (b), as in- 
dicated by the arrows. Panel (c) shows the same trajectories 
as in (b) in the Markov limit x ^ 00. Panel (d) shows the 
symmetry of the most probable trajectories to a state (g,p) 
and from the state (g, —p) in Markov systems in thermal equi- 
librium. 

For a nonparabolic potential U {q) one can find the tra- 
jectory to a given state [q,p) numerically by integrating 
Eq. (24) backward in time from i = 0, using the values 



of 9(07= q and q{Q) — p. To find the optimal escape tra- 
jectory, on the other hand, one can differentiate Eq. (24) 
over time, which leads to equation 

q + xq + {U\[{q) + Chx\q - xU'h{q) = (30) 

One can then seek the solution of this equation by the 
shooting method, starting from t — > —00, where it has 
the form of a sum of exponentials, see Eq. (29). The 
coefficients 51,2,3 have to be found from the condition 
that the trajectory arrives at {qs,P — 0) for t ^ 00 and 
that the solution satisfies the initial intego-differential 
equation (24), which leads to the relation between (71,2.3 
given by the second equation in (29). 



VI. SHORT NOISE CORRELATION TIME 

An explicit solution for the probability distribution 
and the escape rate can be obtained in the case of short 
correlation time of the noise and short delay time of 
the bath, when functions (f>f{t) and 4'h{t) are close to 
(5-functions. To do this it is convenient to eliminate f{t) 
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from the functional (14 1. Then the variational problem 



for reaching a given state {q,p) takes the form 



R{q,p) =min 



dt 



dt'x{t)<j>f{t~t')x{t') 

) 

dtx(t) m + K{m) 



— oo 

t 



+(3 dt'4>h{t~t')q{t') 



(31) 



The variational problem for the activation energy of es- 
cape is given by Eq. (31 ) in which the integrals over t go 
from —oo to oo rather than from — oo to 0. 

If the dissipation has no delay and the noise is white, 
(j)f[t) = ATD5(t) and (l)h{t) = 4TkBT6{t), the optimal 
trajectories for the variational problem ([3l|) are 



[/; - 2rq - 0, 



(32) 



In this approximation, which is of the zeroth-order in 
delay and correlation. 



R^°Hq,p) = ^p^ + Uhiq)-Uhi<la), 
Rf = Uhiqs) - Uhiqa). 



(33) 



If the correlation time of 4'f,(t)h is short compared to 
the relaxation time of the system, one can write 

<j)fMit) = 27r<Pf^h{0)Sit) + [0/,ft(t) - 2^$/,;,(0)5(t)] 

and then consider the term in the brackets as a pertur- 
bation. We choose T and D in the same way as for S- 
correlated (j)f^h{t), 

/oo 
dt(j>f{t) = 4ri:>, 
-oo 

/oo 
dtcj>h{t) ^ ATkBT. 
-oo 

To the first order in the delay/correlation perturbation, 
the correction to R{q,p) can be calculated using the ze- 
roth order trajectories ( [32| ). In the integrals over t' in 
Eq. ( [3I] ) it is convenient to expand x(t') sa x{t) + {t' — 
t)x(t)and q{t') « q{t) + (f - t)q{t). Substituting this 
into Eq. (l3l|) one obtains 



i?(g,p)«i?W(9,P) + r(t/-4)p^ 
tjj;= I t<j)f^h{t)dt/ I (f>f,h{t)dt. (34) 



f t<j)f^h{t)dt / <j}f,h{t)dt. 
JO ' Jo 

The parameters tf and characterize the widths 
of the correlators (j)f{t) and (t>h{t), respectively. From 
Eq. (34), the correction to R{q,p) is generally of the first 
order in these widths (cf. Ref . [5T] where there was found 
a correction due to noise correlations for an overdamped 
system with no delay). The broadening of the noise 
correlator 4>f{t) s-nd the dissipation delay [the broad- 
ening of 4'h{'t)\ act in opposite directions: the larger is 



tf the smaller is the probability (13) to reach a remote 
state, whereas the increase of th increases this probabil- 
ity. In thermal equilibrium the two effects compensate 
each other. 

From Eq. (34), there is no correction to the activation 



energy of escape Ra of the first order in the width of the 
correlators 4)f,h{t), because p = at the saddle point. 
The lowest-order correction appears in the second order. 
It can be calculated similarly to the above procedure, 
taking into account that the integrals over t in Eq. (31) 
now run to 00. In the integrals over t' one should expand 
x(t') and q{t') about x(0 and q{t), respectively, to the 
second order mt — t' . Then Eq. (32) gives 



Ra ~ R'^a 



I' f h — 



t^^f^h{t)dt. 



— 00 
00 



^f,k{t)dt, (35) 



where q(t) is given by Eq. (32 ) with boundary conditions 
q{t) — > qa for t —00 and q{t) — )■ qg for t 00. Again, 
the changes of the escape activation energy Ra due noise 
correlations and dissipation delay have opposite signs and 
compensate each other for thermal noise. 



VII. CONCLUSIONS 

In this paper we have considered large rare fluctua- 
tions induced by Gaussian noise in systems with delayed 
dissipation. The dissipation comes from coupling to a 
thermal bath; the corresponding friction force depends 
on the history of the system motion and is described by 
an integral of the velocity with an appropriate kernel. 
The noise, along with the part that comes from the ther- 
mal bath, can have another source. 

The proposed formulation reduces the problems of 
finding the logarithm of the probability distribution over 
the phase space of the system and the effective activation 
energy of escape from a metastable state to variational 
problems. The extreme trajectories of the respective vari- 
ational functionals provide the most probable paths fol- 
lowed by the system in a fluctuation to a given state or 
in escape and also the most probable corresponding re- 
alizations of the noise. 

We show that, if the noise is coming from the ther- 
mal bath responsible for the dissipation, the logarithm 
of the probability distribution has a familiar Boltzmann 
form, with a renormalized potential due to the coupling 
to the bath. Closed-form equations for the most prob- 
able system trajectories are obtained in this case. We 
describe both the portion of the trajectory leading to a 
state and followed after the state had been reached. It is 
shown that, if the state corresponds to zero momentum, 
these portions are time-reversal symmetric. In contrast 
to Markov systems, however, the phase trajectories are 
smooth at the observation point. This is a clear signa- 
ture of the dissipation delay. Also, in contrast to Markov 
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systems, one may not compare the most probable tra- 
jectory to a given point {q,p) in the phase space of the 
system with the trajectory from this point, since the lat- 
ter is not defined for a system with delayed dissipation. 
Explicit solutions are obtained for the exponentially cor- 
related in time noise and are used to illustrate the general 
properties of the trajectories. 

The case of the noise and dissipation with short corre- 
lation times is analyzed. It is shown that the logarithm 
of the probability distribution has corrections of the first 
order in the widths of the time correlation function of the 
noise and the dissipation kernel. In contrast, the effec- 
tive activation energy of switching has only second-order 
corrections in these widths. 

Delayed dissipation is commonly encountered in physi- 
cal systems, the simplest example being dissipation from 
the coupling to an electromagnetic environment with 



frequency-dependent impedance. This makes it possible 
to test the predictions of the present paper in experiments 
similar to those used to observe optimal fluctuational 
paths in effectively Markovian systems, in particular in 
nano- and micromechanical resonators. Underdamped 
resonators with strong nonlinearity would be advanta- 
geous for such experiments, as there the vibration fre- 
quency strongly varies as a result of energy fluctuations. 
Therefore the density of states of the excitations of the 
thermal bath that resonate with the system and cause 
energy absorption varies as well. Another candidate sys- 
tems of significant current interest are optomechanical 
systems in high-Q optical cavities, where the density of 
states of the electromagnetic bath has strong dispersion. 

The research of M.I.D. was supported in part by NSF 
Grant No. CMMI-0900666 and DARPA through the DE- 
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